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The decay of a crystalline cone below the roughening transition is stud- 
ied. We consider local mass transport through surface diffusion, focusing on 
the two cases of diffusion limited and attachment-detachment limited step 
kinetics. In both cases, we describe the decay kinetics in terms of step flow 
models. Numerical simulations of the models indicate that in the attachment- 
detachment limited case the system undergoes a step bunching instability if 
the repulsive interactions between steps are weak. Such an instability does not 
occur in the diffusion limited case. In stable cases the height profile, h(r,t), 
is flat at radii r < R(t) ~ t l / A . Outside this flat region the height profile 
obeys the scaling scenario dh/dr = ^"(rt -1 / 4 ). A scaling ansatz for the time- 
dependent profile of the cone yields analytical values for the scaling exponents 
and a differential equation for the scaling function. In the long time limit this 
equation provides an exact description of the discrete step dynamics. It ad- 
mits a family of solutions and the mechanism responsible for the selection 
of a unique scaling function is discussed in detail. Finally we generalize the 
model and consider permeable steps by allowing direct adatom hops between 
neighboring terraces. We argue that step permeability does not change the 
scaling behavior of the system, and its only effect is a renormalization of some 
of the parameters. 
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I. INTRODUCTION 

The properties of crystalline nanostructures are of considerable interest, because of the 
technological importance of nanostructures in fabrication of electronic devices. Kinetic prop- 
erties of nanostructures attracted particular attention, since in many cases nanostructures 
are thermodynamically unstable and tend to decay with time. Such decay processes have 
been studied both theoretically and experimentally. 

Considerable effort was devoted to the study of periodic structures. The decay of one and 
two dimensional gratings was studied extensively. The emerging experimental pictureETEl is 
that below the roughening temperature these structures decay in a shape preserving manner. 
Macroscopic facets are observed at the maxima and minima of the gratings. Although these 
systems are out of equilibrium, the appearance of facets is a manifestation of the cusp 
singularity of the surface free energy at the high symmetry crystalline orientation. 

There are basically two theoretical approaches to the problem of surface evolution in 
general and nanostructure decay in particular. On the one hand, tfrfpje are phenomenological 
models which treat the crystal surface as a continuous mediumcTEJ. The evolution of the 
surface is then driven by the tendency of the system to lower its free energy (given in terms 
of continuous spatial variables). The advantage of such models is that they are relatively 
simple and can sometimes lead to analytical predictions of surface evolution. However, these 
models ignore the discrete nature of surface steps, which may become important below the 
roughening transition. In addition, most of them rely on assumptions of small surface slope 
and/or surface curvature (with the exception of Ref. [12]), and are unable to properly treat 
the behavior of the macroscopic facets observed experimentally. 

On the other hand, there are models, which treat surface evolution on a smaller scale. 
Among these are microscopic-mqdels&cfJ, where the basic degrees of freedom are individual 
atoms, and step flow modelstfrM. These models are usually solved numerically, and provide 
results which can be directly related to the microscopic dynamics. However, in most cases, 
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it is difficult to understand the behavior of the system on larger length scales on the basis 
of these results. 

Research efforts were also directed towards isolated surface structures. The de- 
cay of isolatecLrStep bunches, islands or hills was studied both experimentallyQ'B and 
theoreticallytJc3Ej. Here too there is an apparent gap between the microscopic and macro- 
scopic theoretical approaches. 

In this work we attempt to bridge over this gap in the case of a simple surface structure, 
i.e. an infinite crystalline cone. We give a complete account of the surface dynamics based 
on a step flow model, and then derive a continuum model which gives a very accurate 
description of the evolution of the cone and becomes exact in the long time limit. We do 
not make any assumptions of small slope or small surface curvature. 

The crystalline cone consists of an infinite number of circular concentric steps. A similar 
system was studied by Rettori and Villain^ who considered the decay of bidirectional surface 
modulations. Their results are relevant in the case of small amplitude modulations when the 
profile peaks and valleys affect each other. Our work addresses the opposite situation when 
a single peak can be considered as an isolated structure and in this sense is complimentary 
to theirs. 

Below the roughening transition atomic steps have a finite free energy. Their existence 
on the surface strongly affects its morphological evolution. In many cases, one can ignore 
the formation of islands and voids on the surface and consider only adatom diffusion and 
attachment and detachment processes to and from step edges. The decay of a nanostructure 
is then dominated by the motion of steps. In order to describe the decay process mathe- 
matically, one has to solve the diffusion equation for adatoms on the terraces between the 
steps, .with boundary conditions at the step edges, in the spirit of the Burton-Cabrera- Frank 
models. If the geometry of the nanostructure is simple, this procedure leads to a set of cou- 
pled equations of motion for the steps. Our goal in the present work is to construct and 
solve these equations of motion for the simple case of an infinite crystalline cone. A partial 
account of this work is found in Ref. |[ The kinetic step model for the cone is derived in 
section II. In section III we carry out numerical simulations of the model and examine the 
evolution of surface morphology under various conditions. 

Previous experimental and theoretical research on decay of nanostructures has demon- 
strated that in most cases the surface reaches a scaling state where the typical length scale 
depends on time algebraically. Our simulations show that the cone profile also exhibits such 
a scaling behavior. In section IV we show analytically that the step flow model admits 
such solutions. We calculate the scaling exponents and derive a continuum equation for the 
scaling function. The properties of the scaling function are analyzed in sections V and VI. 

It has been suggested that in some materials, steps may be permeable to flow of atoms; i.e. 
atoms can hop between neighboring terraces without attaching to the step which separates 
them. In section VII we discuss the consequences of step permeability on the decay of an 
infinite cone, and show that permeability does not change the qualitative behavior of the 
system. 



II. STEP FLOW MODEL OF A CRYSTALLINE CONE 

We now consider a crystalline surface which consists of flat terraces, parallel to a high 
symmetry plane of the crystal and separated by atomic steps. Islands and voids are ignored- 
The evolution of such systems of steps was treated long ago by Burton, Cabrera and Franktl 
(BCF). The BCF theory assumes that mass transfer between terraces is governed by diffusion 
of adatoms on the terraces. These atoms are emitted and absorbed at step edges, which 
according to BCF act as perfect sinks and sources. This last assumption, however, is valid 
only when adatom diffusion is the rate limiting process. Modifications of the BCF model 
account for the finite rate of adatom attachment-detachment processes at step edges. The 
BCF model was also generalized to include elastic and entropic interactions between steps 



see for example Refs. In this section we construct such a generalized BCF model 



for the morphological evolution of a conic hill on a crystal surface. 
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Consider the surface of an infinite crystalline cone, which consists of circular concentric 
steps of radii rj(i), separated by flat terraces (Fig. [I]). The index i grows in the direction 
away from the center of the cone. 




FIG. 1. Illustration of the step configuration in a crystalline cone. The steps height is exagger- 
ated. 



These steps may absorb or emit atoms, which then diffuse across the neighboring terraces 
with a diffusion constant D s . Assuming no deposition of new material, no evaporation and 
no transport through the bulk, the adatom concentration Cj(r) on the zth terrace satisfies 
the diffusion equation 



DN 2 CAr) 



dt 



where f is a two dimensional vector parallel to the high symmetry terraces. In most situa- 
tions, the time scale associated with step motion is much larger than the time scale of surface 
diffusion. One can therefore assume that the adatom diffusion field is always in a steady 
state; i.e., for any step configuration, the diffusion field reaches a steady state before the 
steps move significantly. Within this quasi-static approximation, the r.h.s. of the diffusion 
equation can be neglected. Using the radial symmetry of the cone, we can write the static 
diffusion equation in polar coordinates as: 



dr 2 

The general solution of this equation is 

CAr) - 



d 2 Q(r) | ldd 
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The coefficients A^ and Bi are determined by the boundary conditions at the step edges. To 
define these conditions, we assume that the flux of atoms at the two step edges bounding the 
zth terrace is determined by first order kinetics, characterized by an attachment-detachment 
rate coefficient k: 



Or 
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where C^ q is the equilibrium concentration of atoms on the terrace adjacent to the ith step. 
Using these boundary conditions, we find that the constants Ai are given by 

* = : . ■ '*), , , ■ W 



In r, - In r i+ i - % f- + — 
Employing mass conservation at the step we obtain the step velocity: 



drj „ ^ / 9Ci <9C } - 
QD 



dt s \ dr dr 



no A — , (5) 



where f2 is the atomic area of the solid. 

In order to complete the solution of the diffusion problem we have to determine the 
equilibrium adatom concentration at step edges. This concentration depends on the local 
step curvature and on the radii of neighboring steps. According to the Gibbs-Thomson 
relation, C^ 9 is given by 

= ^exp|^(l + |) , (6) 

in units where the Boltzmann constant is equal to 1. Here /ij is the chemical potential 
associated with the addition of an atom to the solid at the ith step, T is the temperature 
and C eq is the adatom equilibrium concentration at the edge of a straight isolated step. To 
evaluate the step chemical potential, we take into account the step line tension, T, and 
a repulsive interaction between nearest neighbor steps. The magnitude of this interaction 
is inversely proportional to the square of the distance between tha^teps at large distances. 
Such a dependence is consistent with entropic as well as elasticE3E2l interactions between 



straight steps. We follow Ref. 16 and take the interaction energy between step i + 1 and an 



atomic segment of the ith step to be 



U (r i: r l+1 ) = ■ — -2 , (7) 

(n + r i+1 ) {r i+1 - n) 

where G is the interaction strength. The chemical potential of the zth step is then given by 

This equation applies also to the chemical potential of the first step if we set r = 0. 

As we show below, in the long time limit the distance between steps is small compared 
with the step radius. In this limit we can approximate the step chemical potential by 

, i= ™ +S1G (^l 1 — s ). (9) 

This approximation simplifies the algebra considerably, and we verified that it does not 
affect the analytical and numerical results described below. 

We are now ready to write down the step equations of motion using Eqs. (|), (P) and 
the Gibbs-Thomson relation (^). It is convenient to use dimensionless radii, pi, and dimen- 
sionless time r: 



T 



T \ 2 (_ D R T X -J 



In terms of these variables, the equations of motion take the form 



dpi Clj Oj_X . , 

Pi = — = , with (10) 

dr pi 

6 — 6+1 



(l-g)ln-a--gfi + ^- 
1 / 2p l+ i 1 2pi_! 



Pi \Pi+i + Pi (pi+i-pi) Pi + Pi-i (pi-pi-i) J 

Eqs. ( [TO] ) depend on two parameters: g = measures the strength of step-step inter- 
actions G relative to the line tension T, while the parameter q = (1 + |p^) _1 determines 

the rate limiting process in the system. When q — > 1, diffusion across terraces is fast and 
the rate limiting process is attachment and detachment of adatoms to and from steps. On 
the other hand, when q — > 0, the steps act as perfect sinks and the rate limiting process is 
diffusion across terraces. 



III. RESULTS OF SIMULATIONS 

We integrated Eqs. (|TTJ) numerically both in the case of diffusion limited kinetics (DL) 
and in the case of attachment-detachment limited kinetics (ADL). The initial configuration 
was a uniform step train. In principle, the initial step separation is a parameter of the model. 
However, for any value of this parameter one can change the units of length and time and 
get the same equations of motion with initial step separation of unity and different values of 
the parameters g and q. Thus it is sufficient to consider an initial step separation of unity. 

When the repulsive interactions between steps are weak (i.e. g is small), there is a striking 
difference between the dynamics in the DL and ADL limits. In the ADL case the system 
becomes unstable towards step bunching, whereas in the DL case there is no such instability. 
However, when g is large enough the instability disappears even in the ADL case. Let us 
first discuss situations where the step bunching instability does not occur. Fig. ^| shows the 
time evolution of the system in the ADL and the DL cases with a relatively large value of g. 
Each line in the figures describes the evolution of the radius of one step. We note that the 
innermost step shrinks while the other steps expand and absorb the atoms emitted by the 
first step. When the innermost step disappears, the next step starts shrinking and so on. 
Our observations indicate that the time at which the nth step disappears, r n , grows with n 
as r n ~ n 4 . This process results in a propagating front, which leaves a growing plateau or 
a facet at the center of the cone. At large times, the (dimensionless) position of this front 
behaves as /Witf^) ~ t 1//4 . This is shown by the dashed lines in Fig. ^|. 

This power law is an indication of a much more general and interesting phenomenon. It 
turns out that at large times, not only the front position but also the positions of minimal 
and maximal step densities scale as r 1//4 . In fact, the step density D(p,r), defined as the 
inverse step separation, obeys the following scaling scenario: There exist scaling exponents 
a, (3 and 7, which define the scaled variables x = pr"^ 7 and 9 = r 1 /" 7 . In terms of these 
variables 

D(j>,T) = e a F(x,9), (11) 

where the scaling function F is a periodic function of 9 with some period 9q. Our ansatz is 
somewhat weaker than standard scaling hypotheses, which would assume the scaling function 
F is independent of 9. We introduce this periodic dependence because our simulations 
strongly indicate a periodic behavior, generated by the motion of the first step (see Fig. ||) . 
Thus the disappearance time of step n is r n = (n#o) 7 . An immediate consequence of the 
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scaling ansatz is that if we define 9 = 9 + h9q with < 9 < 9q, and plot 9 a D (p, r) against 
x, all the data with different values of n and the same value of 9 collapse onto a single curve, 
F(x,6). 




x 10 




x10" 



FIG. 2. Time evolution of the step radii in the (a) ADL and (b) DL cases with g = 0.01. The 
radius of the facet edge can be fitted by a t 1 / 4 power law (dashed lines). 

To verify that our system obeys this scaling ansatz, we define the step density at a 
discrete set of points in the middle of the terraces: 



D 



'Pi{r) +Pi+i(r) 



1 



p m (r) - pi(r) ' 



(12) 



Fig. |3| shows plots of D(p, r) as a function of x = pr^ 1 ^ 4 for a fixed value of 9 and 6 different 
values of n in the ADL and the DL cases. The excellent data collapse shows that our scaling 
ansatz indeed holds with a = 0, (3 = 1 and 7 = 4 in both cases. 




FIG. 3. Data collapse of the density function in the (a) ADL and (b) DL cases with g = 0.01. 
The values of the scaling exponents used here are a = 0, f3 = 1 and 7 = 4. These figures show 
density functions with 6 different values of n and the same value of 8, as a function of x = pr -1 / 4 . 
Different symbols corresponds to different values of n. The unsealed data is shown in the insets. 

We now turn to discuss effects of interactions between steps. As shown above, the 
behavior of the cone in the ADL and the DL limits is very similar when the repulsion 
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between steps is strong. We have already mentioned that when the repulsion is weak, the 
behavior of the system in the DL limit is very different from its behavior in the ADL case. 
This is not surprising, since linear stability analysis (see Appendix A) of the two cases in 
the absence of step-step interactions predicts that the ADL case is unstable towards step 
bunching while the DL case is marginal. In intermediate cases (0 < q < 1) the system is 
unstable when the interactions are sufficiently weak. 

Fig. f| shows the time evolution and the data collapse of the density function in the 
DL case when step-step interactions are weak. As one can see, the behavior of the cone 
is qualitatively similar to the large g example (Fig. |2] (b)) and that the scaling ansatz still 
holds. Quantitatively, the step density near the facet edge is much higher for small values 
of g. Also the dependence of the scaling function on scaled time within a collapse period is 
much more pronounced when the interactions are weak. 



a) b) 




T *io 5 x=px- 1/4 

FIG. 4. (a) Time evolution of the step radii in the DL case with g = 1CP 6 . (b) Data collapse of 
the density function in the same system. Different symbols correspond to different times. 

The dependence of step kinetics on the strength of the interactions is much more com- 
plicated in the ADL case. Fig. |5] shows the evolution and when possible, the scaled density 
function of a series of ADL systems which differ in the value of g. The two quantitative 
observations we made in the DL case hold also here: The high step density region near the 
front becomes more dense, and the time dependence of the scaling function becomes more 
pronounced as the value of g is reduced (Fig. |5] (a)-(d)). In addition, we noted that the 
approach of the system towards the scaling state is slower. Adjacent to the dense region at 
the edge of the facet is another region of very low density of steps. This region becomes less 
dense as the value of g is reduced. Below a critical value, g c , of the interaction parameter, a 
few steps between the low density region and the facet edge form a bunch. At this stage the 
scaling ansatz breaks down. The system no longer exhibits a simple periodic nature, which 
results from the collapse of single steps, one at a time. Instead, it seems to adopt a compli- 
cated almost periodic pattern which involves the collapse of many steps in each period (Fig. 
[5] (e)). Hints for this periodicity changes are already present in Fig. H (a) and (c), where the 
steps crossing the low density region follow a threefold periodicity. lithe value of g is further 
reduced, bunches of steps collapse together, and finally the whole system becomes unstable 
toward step bunching (Fig. [5] (f)). Neighboring steps merge and form a bunch, which in 
turn merges with another bunch and so on. Step bunches rather than isolated steps become 
the dynamic objects in the system. 



7 




b) 



x 10 



2.5 



1.5 



0.5 



d) 



2.5 



1.5 



0.5 




x=px 



■1/4 





x 10 



FIG. 5. Time evolution and scaled density functions of ADL systems with different interaction 
strengths, (a) and (b): g = 10" 3 , (c) and (d): g = 5 • 10" 4 , (e): g = 2.5 • 10" 4 and (f): g = 1(T 6 . 



IV. SCALING ANALYSIS AND THE CONTINUUM MODEL 

The above results suggest that when the step bunching instability does not occur, the 
time evolution of the system can be described by a step density function, which is continuous 
in both position and time variables. In this section we derive such a continuum model by 
carrying out a scaling analysis. We obtain an equation which governs the scaling function 
and evaluate the scaling exponents analytically. 
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Motivated by the simulation results, we assume that at long times the scaling ansatz, 
Eq. ([TIP , holds. This together with conservation of the total volume of the system, already 
determines the two scaling exponents a and f3. First, we derive a relation between these 
scaling exponents by considering the height profile h(p,r). Assuming steps of unit height, 
the height difference between two points on the surface is given by the number of steps 
between them. The continuous analog of this statement can be used to derive the following 
relation between the height profile and the step density: 

h(p,r) = ho(r) - [ P D(p',T)dp' , (13) 
J o 

where ho(r) is the height at the origin. Far enough (in the limit p — > oo), h(p,r) does not 
change with time. Calculating the height change at infinity using Eq. (|13D we find 

hoir) - WO) - r~ / [Fix, 9) - F(oo, 0)1 dx = , (14) 

Jo 

where we have changed the integration variable and used the definition (Eq. (|TTD ) of the 
scaling function F. On the other hand, h (0) — ho(r n ) = n because r n is the time of 
disappearance of the nth step. 7 satisfies the relation r n ~ n 7 , and therefore we have 
/io(0) — ho(r n ) ~ ■ This and the r^ a+/3 ^ 7 dependence in Eq. ([14]) lead to the relation 
a + (3 = 1. 

In addition, conservation of the total volume of the system, V, implies that 

V = 2tt / ph(p,r)dp (15) 
Jo 

is independent of r. Integration by parts of the derivative of this integral with respect to r 
yields the following equation: 







or 



(16) 



Evaluating this integral in terms of the function F and the scaled variables x and 9 we 
obtain the equation 



x 



aF(x,6) dF(x,\ 
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dx = 



(17) 



This can be satisfied for all 9 only if a = 0, since F is a periodic function of 9. Combining 
this result with the previously obtained relation a + (3 = 1, we conclude that 

a = , p = 1 . (18) 

Thus, 7 is the only nontrivial scaling exponent in the model. To evaluate 7 and the 
scaling function F, we continue with the equation for the full time derivative of the step 
density D: 

<W__dD_ dD_ dp 

dr dr dp dr 

Eq. flHSP can be evaluated in the middle of the terrace between two steps (i.e. at p = 
(Pi + Pi+i)/2). The l.h.s. of (|T9| ) is calculated by taking the time derivative of Eq. (|T2"D: 
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dD/dr = —D 2 (p i+ i — p\). This together with the fact that dp/dr = (p\ + p i+ i)/2 leads to 
the relation 



+ -A> = 0. (20) 

Now we change variables to 9 and Xi = pi9~ 1 (since (3=1), and transform Eq. fl2T)| ) into 
an equation for the scaling function F: 

The step velocities pi and can also be expressed in terms of the x^s, but we defer this 
algebraic manipulation to a later stage. 

Our next goal is to take a continuum limit of Eq. ([H]) in the variable x = (a^+i + Xi)/2. 
Such a continuum limit becomes exact in the long time limit. To see this, let us rewrite Eq. 
(12]) in terms of Xj's, 9 and F: 

F((x i+1 +Xi)/2,9) 

According to this equation, the difference between successive x^s is of order 9~ l wherever F 
does not vanish. In the large 9 (long time) limit these differences become vanishingly small. 
This justifies the continuum limit in the variable x. 

In practice, we take the continuum limit in the following way. We evaluate the function 
F at the position (x i+k + x i+k+ i)/2 by using it's Taylor expansion 

Xj+k + Xj+k+i \_ 9- 1 _ ™ 1 d n F(x,9) / x i+k + x i+k+1 \ n 
2 ' J ~ x i+k+1 - x l+k ~ ^ Q n\ dx n \ 2 xj ■ 

As long as k is finite, the difference Xi +k — x is small in the long time limit. It is therefore 
useful to expand x i+k around x 



%i+k 



X 



n=l 



(24) 



and insert this expansion into Eq. (f^J). The resulting equation is expanded as a power series 
in 9~ x . By requiring that the equation be obeyed to all orders in 9~ 1 , we can calculate the 
coefficients (f>k n for any desired k and n. These coefficients will involve the function F and 
its derivatives with respect to x, which are all periodic functions of 9. 

Next, we express the velocities pi and p\ + \ in terms of the scaled radii Xi-2, Xj-i, . . . , 
x i+3 using Eq. (TTD|) and the transformation to scaled variables. We then use Eq. ( [Zip to 
expand the velocities in powers of 9~ x . The result of this expansion is 

p l+1 +p l = 9- 3 (A + 0(9- 1 )) , 

p i+1 -p l = 9-\B + 0{9- 1 )) . (25) 

A and B are known expressions involving F, F', F", F'", F"", where the primes denote 
partial derivatives with respect to x. The full expressions are given in Appendix B. The 
existence of derivatives up to fourth order in this equation is a consequence of the fact that 
each step "interacts" with four other steps (two on each side) through the equations of 
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motion (|T0|). Inserting Eqs. fl25|) into Eq. (0), we obtain the following differential equation 
forF: 



— F' • — + 9^ (F'-4 + F 2 b) + -~ + 0(9^ 5 ) = . (26) 
7 V 2 / 7 39 



- +IP -'(F'4 + F'B) + e - 9F 
7 V 2 / 7 

Consider Eq. (KB]) at large Our expansion in the small parameter # _1 is valid only at 
values of x where F does not diverge or vanish (see above). Therefore, the first term in Eq. 
( p6|) is 0(1). This term has to be canceled by the second term if we require F to satisfy a 
single differential equation. Hence, we must have 

7 = 4. (27) 

The fourth term vanishes as 9 — > oo since 7 — 5 < 0, and the third term must vanish as well. 
Therefore, in the large 9 limit, F is only a function of x, and we are left with an ordinary 
differential equation for F: 

F'(^-|)+F 2 B=0. (28) 

The detailed form of this equation in the DL and ADL cases is also given in Appendix B. 

Let us emphasize several important properties of our scaling analysis. First, the values 
of the scaling exponents we calculated (a = 0, (3 = 1 and 7 = 4) are consistent with the 
results of numerical simulations (see above). Secondly, our continuum model is valid for 
arbitrarily large surface curvature and slope (unlike other treatmentsE-M) . Moreover, since 
our model is an expansion in the truly small parameter 9~ x (see Eq. fl2~2|) ) it becomes exact 
in the large 9 (long time) limit. Finally, note that in going to the continuum limit we lost 
the periodic dependence of F on 9. This periodicity is generated by the first step which 
follows a unique equation of motion. We did not incorporate this unique behavior into the 
continuum model and therefore should not be surprised that this information is lost. 

As we emphasized, our continuum model is an exact representation of the original discrete 
system. For this reason it is interesting to compare it with other continuum models, which do 
not emerge as limits of discrete systems of steps. Many authors use the continuity equation 

§ + VJ = 0, (29) 

to account for the surface evolution. Here h is the surface height and J is the current 
density of diffusing adatoms. This equation is of course correct and reduces the problem to 

calculating J. It is widely assumed@Q that J is proportional to the gradient of the surface 
chemical potential /i, as expected in diffusive systems. In the case of attachment-detachment 



limited systems it was suggested by Nozieres£3 that the chemical potential gradient should 
be divided by the profile slope. Our model is consistent with the first picture in the DL case 
and with the second in the ADL case. 

To show this we take the gradient of Eq. (|29|). This leads to 



dD 

-^ocVV-J, (30) 



where — D is the gradient of the profile. We now return to Eq. (|28|) and transform it back 
to an equation for the step density _D(p, r). Remembering that the term — ^f- in Eq. 
arises from the time derivative of D, we find that in the general case we can write 
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Eq. (|3l|) is nothing but the radial component of Eq. fl3"0|) written in terms of the dimensionless 
variables p and r with 



l-q + 2qD dp V ' 

Therefore in the DL case (q = 0) the adatom current is indeed proportional to the 
chemical potential gradient. In the ADL case (q = 1) the chemical potential gradient is 
divided by the local step density as suggested by Nozieres. In addition we note that in 
the limit p — > oo (where the steps are nearly straight), our chemical potential (Eq. fl3"I|)) 
becomes identical to the chemical potential of Ref. [12]. 



V. PROPERTIES OF THE SCALING FUNCTION 



In this section we use Eq. Q2"8p to study properties of the scaling function F. We also 
discuss the boundary conditions necessary to solve Eq. fl2"8"|). We begin by noting that 
according to the numerical simulations of the step model, there is a growing plateau or facet 
at the top of the cone. Our numerical solutions for the scaling function indicate that all 
physical solutions indeed have such a special point, which can be identified as the facet edge. 
Let us examine the behavior of the step density on the facet and at its edge. At any given 
time there is only one single step on the facet (the first step during its collapse towards the 
origin). We have seen that the size of the facet grows indefinitely, and therefore the step 
density on the facet vanishes in the long time scaling limit. Moreover, in Appendix C we 
show that in the long time limit, the step density is a continuous function even at the facet 
edge; i.e., it goes to zero continuously as the facet edge is approached from above. 

Denoting the scaled position of the facet edge by x , the above observations can be 
expressed in the form 

F(x) = , Vx < x . (33) 

Note that F(x) = is not a solution of Eq. (^). This does not contradict the continuum 
model, since the model was derived only for the case of a finite step density (see Eq. (|2~2"D ). 
Thus the scaling function has to satisfy Eq. (^) only at x > xq, and x = xq is a singular 
point. 

Let us now study the nature of this singularity at the facet edge. Returning to our 
simulation data we note that near the facet edge, the scaling function is extremely steep 



see Fig. 0). Indeed, the expansion of Eq. (28) in small F suggests that in the vicinity of xq, 



F(x) can be written as a power series in ^Jx — Xq. Thus although F is continuous at xq, all 
its derivatives with respect to x diverge at the singular point. 

We now turn to discuss the boundary conditions necessary to solve Eq. (p8|). First, we 
examine the behavior of F at infinity. Far enough from the origin the steps do not move. 
Hence in the limit p — > oo the step density remains at its initial value. This implies that 

lim Fix) = 1 . (34) 

x— >oo 

It turns out that this condition corresponds to two boundary conditions, since it is only 
possible to satisfy Eq. (0) if F' vanishes at infinity. In fact, by considering the asymptotic 
expansion of F it can be shown that for large x 



12 



F = 1 + a 4 x~ 4 + O (x~ 8 ) , (35) 

where — 3 (1 + g) in the DL case and oj = | (1 + g) in the ADL case. 

In order to solve Eq. (|28|) we need two additional boundary conditions at x = Xq. The 
first of these is F(xq) = (see Eq. (|33"D). Another boundary condition at x = xo can be 
derived by enforcing volume conservation. From Eq. (|l"5|) , the volume change during the 
time t is given by 

roo 

2tt p[h(p,r)-h(p,0)}dp. (36) 
Jo 

Integrating this equation by parts and requiring volume conservation, we get 

roo 

/ p 2 [D(p,r)-D(p,0)}dp = 0. (37) 
Jo 

In terms of scaled variables the last equation can be written as 

roo nr3 

x 2 (F - 1) dx = / x 2 (F - 1) dx - = , (38) 



.ro 



where we have used the fact that F vanishes below Xo- The evaluation of the integral in Eq. 
(|38|) can be done by multiplying Eq. ( |28|) by x 2 and integrating it with respect to x. The 
resulting equation is 

x i(Lj^ + F 2B)dx= —dx. (39) 
xq y 2 J J xo 4 

The integral = I x 2 (j-^ + -F 2 £>) <ix on the l.h.s. is carried out in appendix C, and is 

expressed in terms of F and its derivatives. The result of this integration combined with 
integration by parts of the r.h.s. of Eq. leads to 



f°° 1 r 

/ x 2 {F - l)dx = - x 3 (F - 1) - 4M 2 

J XQ 3 ^ 



(40) 

x 



Inserting this relation in Eq. (|38|) we obtain the boundary condition 

M 2 \ x =0. (41) 



We have used the facts that the surface term at infinity in Eq. ( [40]) vanishes and F(xq) = 0. 

At this point we have four boundary conditions, two at infinity and two at xq. We may 
now obtain a unique solution of Eq. pS|) if we know the value of xq- What determines 
Xo? To answer this question, considerthe height of the cone at the origin, ho, at time 
t = r n -\ and at time r = r n . r n _i and r n are the disappearance times of two successive 
steps, and therefore ho (T n _i) — ho (r n ) = 1. We can also calculate this difference from Eq. 
(11). Combining these two results we arrive at the relation 



r l/7 _ r 



iL\)f(F-l)dx = -l. (42) 



Recalling that r n = (6orif and using the fact that F vanishes below xo, we can rewrite the 
last equation as 
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poo 

/ (F — 1) dx — x = —0 

JXQ 



(43) 



The integral in Eq. ([4~3"1) can be evaluated by integrating Eq. fl28f) with respect to x. The 
result is: 

poo 

/ (F - 1) dx = AM \ xo + x , (44) 

where the integral Ai = J (j^r + F 2 B s j dx is carried out in Appendix C, and is expressed 
in terms of F and its derivatives. Combining Eq. ( f£3"D and Eq. (fP|) we obtain the relation 

4-MoL = • (45) 

The l.h.s. of this equation depends on x , thus relating x to 9 , the scaled collapse time 
period of the steps. What is then the value of 8q ? It is determined by the motion and 
collapse of the first step. The behavior of the first step is different from that of all the other 
steps, since it does not have neighboring steps with smaller radii. Thus, our continuum 
model, which treats all the steps on equal footing, does not contain any information on 
the value of 8q. We therefore expect a family of valid scaling functions consistent with the 
continuum model with different values of Xq or 6 . The unique value of Xq observed in 
simulations is determined by the discrete nature of the steps, and at this stage we are not 
able to calculate it. 



VI. NUMERICAL EVALUATION OF THE SCALING FUNCTION 

We now find the scaling function numerically. Consider first the DL case. We choose a 



value of x and solve equation (|28| ) starting at x = x . As explained in the previous section, 



F can be expanded in powers of yjx — xq in the vicinity of xq\ 

CO 

F(x) = J2 k n (Vx~^ Y , (46) 



n=l 



where we have already used the boundary condition (p3|). We therefore have three additional 
free parameters in this expansion, which should be determined by the boundary conditions 
(34) and (|41"|). We choose these parameters as the first three odd coefficients in the expansion 
(46): k\, &3 and k§ (the first two even coefficients vanish). Using the boundary condition 
(fl), we can express k 5 as a function of k\ and k% through the relation 



k\ k\ k 3 



k 5 = — - - — + 

6xq 2ki 9x Qgx'lki 

thus reducing the number of free parameters to two. 

Finally we use the boundary condition ( |34T ) to find the coefficients k\ and k^. Since ( |34"D 
addresses the value of F at infinity, it cannotbe easily applied to the expansion of F in the 
vicinity of Xq. We therefore start from an initial guess for the values of k\ and k 3 , solve Eq. 
(f28|) numerically for this choice of parameters and then tune k\ and k^ until the boundary 
condition ( p4|) is satisfied. For a given choice of k\ and k 3 , we found the solution of Eq. 
( P%D in the following way. Numerical integration starting at xo is impossible because the 
derivatives of F diverge there. Therefore, we first used the expansion (f4"ED, truncated at a 
high enough order, to evaluate the function F and its first three derivatives at x = xq + Sx, 
for some choice of 5x. Then, using these values we integrated Eq. (|28|) numerically from 
x + 5x. We made sure that the solution is not sensitive to the choice of 5x. 
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FIG. 6. Numerical solutions for the DL scaling function compared with simulation data. Results 
for three different values of the step-step interaction parameter are shown: (a) g = 0.1, (b) g = 0.01 
and (c) g = 0.001. When the step-step interaction is strong, the Xq solution (dashed line) agrees 
very well with the simulation data (circles). As the value of g is reduced, the best fit solution of 
Eq. ( |28D (solid) deviates from the Xq solution to higher values of xq. 

The above procedure was employed to generate the scaling function for a range of values 
of the scaled facet edge position, Xq, and the following picture emerged. There is a special 
minimal value of Xq, which we denote by Xq. For x < x o there is no solution of Eq. Q28D 
which satisfies the boundary conditions. For each value of Xq > x* Q there is a unique solution 
which satisfies the boundary conditions. Thus, as we anticipated, there is a one dimensional 
family of scaling functions parameterized by the value of xq > Xq. 

Since xq is related to 0q through Eq. (KD, this family of solutions can alternatively 
be labeled by the value of 9q. We used Eq. [E5|) to calculate 9q{xq) and found that it is a 
monotonically decreasing function. In particular, #0(^0) is maximal at Xq. This is reasonable, 
since it means that a smaller facet corresponds to a longer period and thus to a slower facet 
edge. 

Despite the existence of many possible scaling functions, our simulations suggest that 
the system reaches a unique scaling solution independent of initial conditions. What is the 
selected solution? In Figs. M we compare the calculated scaling functions with simulation 
data in the DL case. We do this for three different values of g, the interaction strength 
parameter. When g is large there is an impressive agreement between the simulations data 
and the Xq solution (Fig ^] (a)). When g is reduced (i.e., for weaker interactions between 
steps) the observed scaling function deviates from the Xq solution (Fig. ^| (b) and (c)). 
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However, in these cases there is another solution with a larger value of xq that best fits the 
simulation data. The agreement between this best fit solution and the simulation data is 
again excellent. 

The above observations suggest that Xq is the selected solution in the large g limit. We 
propose the following argument to support this scenario. Since the parameter g is a measure 
of the strength of the step-step interaction G relative to the step line tension T, the large g 
limit is equivalent to the small T limit. The collapse of the first step is driven by the step 
line tension. In the large g limit the collapse driving force is minimal and the collapse time 
period is maximal. As we mentioned above, a long collapse period is equivalent to a large 
value of #0) which corresponds to a small value of Xq. Thus the large g limit corresponds to 
the minimal value of xo, i.e. x = Xq. 



a) b) 

I 1 1 1 1 1 1.5 




FIG. 7. Numerical solutions for the ADL scaling function compared with simulation data. Re- 
sults for two different values of the step-step interaction parameter are shown: (a) g = 0.2 and (b) 
g = 0.01. When the step-step interaction is strong the simulation data (circles) agrees with the 
x% = 1.09 ± 0.07 solution (dashed line). As g is reduced, we notice that the best fit solution of Eq. 
( ]2q ) (solid line) deviates from the Xq = 1.4 ± 0.07 solution to higher values of xq. 

Now consider the ADL case. Although in this limit we can also expand F in yjx — Xq, the 
numerical procedure described above is not an effective method of solution in this case. It 
turns out that the resulting scaling functions are sensitive to the choice of 8x. We therefore 
had to use a different method to solve Eq. fl28|) in the ADL case. Let us denote by x pea k the 
minimal value of x, which corresponds to a local maximum of F . Such a point must exist for 
any continuous function satisfying the integral condition fl3§D and the boundary conditions 
(j33|) and fl34|). By tuning the values of F(x pea k), F"(x pea k) and F"'(x pea k) we found a one 
dimensional family of scaling functions, which satisfy the boundary conditions. This family 
can be parameterized by the value of x pea k or alternatively by the value of the resulting xq . 
As in the DL case, we find that there is a minimal value of Xq, denoted by Xq, below which 
there are no solutions satisfying the boundary conditions. 

Figs. ^ show the simulation data compared with calculated scaling functions for two 
values of the interaction strength g. We again find that for large g the simulation data 
agrees with the Xq solution (dashed line). When g is reduced, the agreement deteriorates 
and there is a different solution with xq > x$ (solid line) which best fits the simulation data. 

VII. EFFECTS OF STEP PERMEABILITY 

The step flow model we introduced in section II assumes that steps are impermeable; 
i.e., adatoms cannot hop between neighboring terraces without being incorporated in the 
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step separating them. In general, however, steps may be permeable. For example, in a 
recent papeiliS, Tanaka et al. interpret their experimental results as evidence that steps 
on Si(OOl) are permeable. We therefore ask the following question: What is the effect of 
step permeability on decay of nanostructures in general, and on the decay of an infinite 
cone in particular? In this section we show that the scaling exponents of an infinite cone of 
permeable steps are identical to the exponents associated with impermeable steps. Moreover, 
the only effect of step permeability on the differential equation for the scaling function is a 
renormalization of one of its parameters. 

Following Ref. [16| we generalize our step flow model to include step permeability. We 
assume that flux ofadatoms between two neighboring terraces due to direct hops is deter- 
mined by first order kinetics. Introducing the permeability coefficient p, we rewrite Eq. @ 
as 

= k(c t \ n -c:«)+p(c l \ n -c l ^Q 

ri 

= k{C^\ n -Ct q )-p{C t \ n -C t ^\ r ) . (47) 

n 

Step permeability does not affect the diffusion equation itself, and therefore the general form 

of the diffusion field, given by Eq. (|2|), remains unaltered. However the coefficients A4 and 

Bi are affected, and the equations for coefficients of different terraces become coupled. As a 

result, the step equations of motion cannot be written in closed form. The scaling analysis 

therefore becomes more cumbersome. 

As in the case of impermeable steps we assume that in the long time limit the scaling 

ansatz Eq. (|i~T| ) holds, with the same definitions of the scaled variables, x and 9 in terms of 

p and r. We slightly change the definition of the dimensionless time r to r = t/to, where to 

is a time scale which we can choose to our convenience. The derivation of the exponents a 

and (3 in section IV did not involve the dynamics of the system. Thus the same derivation 

still holds and the values of a and (3 are not affected by step permeability. We proceed to 

evaluate the exponent 7 and the scaling function F by studying Eq. ([47D . 

For convenience, we parametrize the solutions of the diffusion equation by concentration 

differences instead of the parameters Ai and Bi. These concentration differences, Ui and Vi, 

are defined as follows: 

Ut = C t (r t ) - C e « = A t lnr, + B { - 

Vi = Ci(n) - Ci(r i+1 ) = Ai In — . (48) 

r i+l 

In the continuum limit, Ui and Vi are continuous functions of r and t or alternatively of x 
and 9. We associate the values of Ui and Vi with the middle of the ith terrace, namely 

U{^^,9) =U t , 

v(^±p^,e)=v t . (49) 

The scaling scenario for the functions U and V is 

u{x,9) =e»(u(x) + o(e- 1 )) , 

V (x, 9) = 9 V (y (x) + O (V 1 )) . (50) 
In terms of U, Vi and the dimensionless radii pi, Eq. ([47D takes the form: 



D, 



da 



Or 



-D, 



Or 
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TDM 



i- z = k{u i - C e %) +p(U i - {Ui-! - Vi- X )) 



Pi+i 



QTpi In 

Evaluating Eq. (ETJ) at x — (xi + x i+ i)/2 we can now employ Eq. 
Pi In and pi In in powers of 9" 1 : 



(51) 

to expand £ i5 



Pi 



& = g +9 (£ + s™<)) s - + (r") 



Pi In 



Pi 



Pi+1 

Pi in^ = -i + o(e~i) . 



(52) 



Pi F 

The Xj's dependence of U and V is unknown so we expand: 



u 



Xi—i -\- Xi 



Xi—1 ~t~ Xi 



E 

n=0 

oo 

E 

n=0 



v 1 d n u(x) (Xi_x + x 



n! <9x n 



1 9 n v(a;) / + x 



x 



n\ dx r 



(53) 



Since the difference between successive x's is of order Eq. (|53|) is also an expansion in 
this small parameter. 

Using Eq. (|53|) together with Eqs. Q49|) and (|50|) we isolate v(x) in Eqs. (|5l|), keeping 
only the lowest orders in 9~ l . 



9 1+u v(x) 
9 1+u v{x) 



nrk (e 1+ »u(x) - c eq (l + g(^ + 3FF' 

D S TF + pttT 
QTk (O l+ »u{x) - C eq (~ + g(^ + 3FF' 
D S TF + (p + k)ttT 



(54) 



Since v(x) cannot vanish identically, the two expressions above are consistent only if u < — 1, 
p = — 1 and 



u(x) = C eq 



- + g(— + 3Ff) 

X \ X / 



(55) 



Thus, both sides of Eqs. (|54j ) decay in time. We now subtract the second line of Eq. (|5T 
from the first line and obtain the following equation: 



D S T ( Vi-i V t 



fir \ Pi In ^ ' ftln^ 

V Pi r Pi+1 

Again we expand to lowest order in 6* _1 and isolate v(x): 



(k + 2p) {Ui-Ui-x + V^x) 



(56) 
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(57) 



v(x) does not depend on 9. This implies that v = —2 and v(x) is proportional to u'(x). 

To finish the scaling analysis we return to Eq. (pi)) and find the leading orders in 
Pi + i + pi and — pi. Using Eqs. @ and (]48| ) we find that 



of 



D S T 2 

nr 2 



Vi 



Vi-i 



pi In pi In 



(58) 



Putting together Eqs. ([50]) ([52] ) and fl5"8|) we conclude that to lowest order in 9' 

p l+1 + p l = 6-\A p + 0{6- 1 )), 

p l+l ~p l = e-\B p + o{e- 1 )) . 



(59) 



Ap and Bp are expressions involving F 



F', F", F'" and F"" . Note that the orders of 9' 1 
in these expressions are identical to those in the equivalent expressions in section IV (Eq. 
( p5]) ). These orders are responsible for the result 7 = 4 in the case of impermeable steps. 
Therefore, in the permeable case we also have 7 = 4. 

After expanding p i+1 + p\ and p i+ i — p\ we can use Eq. to obtain the differential 
equation for the scaling function in the permeable case. However this exercise is not nec- 
essary. We avoid it by arguing that the resulting differential equation in the permeable 
case is equivalent to the one in the impermeable case with renormalization of some of the 
parameters. To see this we note that the treatment we presented here is valid also in the 
impermeable case. Thus when p = 0, the general differential equation must be equivalent to 
the equation derived in Section IV. In addition, the attachment/detachment rate k and the 
permeability p affect the step velocities only through the function v(x), which depends on k 
and p only trough the sum k + 2p. Therefore, k and p affect the differential equation itself 
only trough the sum k + 2p. We conclude that the scaling function in the case of permeable 
steps is identical to the scaling function of impermeable steps with k + 2p in the former 
replacing k in the latter. 



VIII. SUMMARY AND DISCUSSION 

We have studied the relaxation of an infinite crystalline cone below the roughening 
temperature, in terms of a step flow model. The model was solved numerically and two types 
of dynamical evolutions were found. When the repulsive interactions are strong enough, 
the decay of the cone proceeds through the collapse of the innermost steps, one step at a 
time. Weak interactions lead to a step bunching instability (except in pure diffusion limited 
kinetics), and the decay process becomes much more complicated and involves collapse of 
bunches of steps. 

Focusing on stable cases, we found that in the long time limit the decaying step system 
obeys a scaling scenario. The step density (i.e. the slope of the height profile), defined as the 

inverse step separation, scales in time according to D(p,r) = t 

function of the scaled position x = r _ ^ 7 p and exhibits a periodic dependence on the scaled 
time 9 = t 1 / 7 . In particular, the position of the facet edge at the top of the cone grows as 
T Ph _ The values of the scaling exponents which fit our simulations are a — 0, f3 — 1 and 
7 = 4. 

Following this observation we used a scaling ansatz to transform the discrete step flow 
model into a continuum description of surface evolution. The basic predictions of this contin- 
uum model are the values of the scaling exponents (which agree with the simulation results) 
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and a differential equation for the scaling function F. This continuum model becomes exact 
in the long time scaling limit, and it breaks down whenever the step density vanishes. This 
fact can be seen both from the derivation of the continuum model and from the resulting 
differential equation, which predicts a singular behavior at the zeros of the scaling function 
F. 

We showed that each physical solution of the continuum equation has a special point, 
x , at which the scaling function, F, vanishes. F is singular at x = x , and we identify this 
special point as the scaled position of a facet edge. Thus the edge of a macroscopic facet in 
the discrete system is a singular point of the scaling function F . 

A detailed analysis of the discrete system revealed a sufficient number of boundary 
conditions, which define a unique scaling function, F, for a given value of the scaled facet 
edge position, xq. However, we were not able to find a unique value for xq, and were left 
with a one dimensional family of solutions, parameterized by xq. A numerical solution of 
the differential equation for F confirmed the existence of this family. We found that there is 
a minimal value of xq, which we denoted by Xq. For every xq > Xq there is a unique scaling 
function, while for Xo < Xq there are no solutions, which satisfy the boundary conditions. 

A comparison of the numerical solutions of F with results of simulations of the discrete 
system leads to the following picture. When the step-step repulsion is strong, there is a 
remarkable agreement between the xq = Xq solution and the simulation data. When the 
magnitude of the step-step interaction is reduced, the system reaches scaling solutions with 
xq > Xq. We therefore advance the hypothesis that in the strong interaction limit the system 
approaches the minimal x^ solution. 

Although our work provides a detailed account of the decay process of an infinite cone, it 
leaves a few unresolved issues. First and most important is the fact that we were not able to 
uniquely determine the value of xq. We did show that xq depends crucially on the detailed 
collapse process of the innermost step, which is not included in the continuum model. Thus, 
in order to evaluate xq, one has to deal directly with the discrete step system. Another open 
issue related to the behavior of the innermost step is the periodic behavior of the scaling 
function. This periodicity is observed in the kinetics of the discrete system of steps, and is 
absent from the solution of the continuum equation. Lastly, when the repulsive interaction 
between steps is weaker than a certain threshold, the system becomes unstable and step 
bunches are formed. The critical value of the interaction below which this instability occurs 
and its dependence on the kinetic parameter q have not been studied so far. We intend to 
address these open questions in future work. 

Finally, we remark that the scaling behavior of the cone profile, predicted in this work, 
is robust. In particular, the existence of a facet growing as r 1//4 , as well the existence of 
a scaling state does not depend on the detailed form of the repulsive interactions between 
steps. A quantitative change of these interactions may alter the scaling function, F, but 
the scaling exponents do not change. Another manifestation of the robustness of the scaling 
solution is the effect of step permeability. In principle, step permeability could have changed 
the scaling behavior of the system entirely. However, we showed that its only effect on the 
scaling solution is to modify the attachment-detachment rate coefficient, k, to k + 2p. 

We are greatful to N. Bartelt, H. C. Jeong, D. J. Liu and J. D. Weeks for helpful 
discussions. This research was supported by grant No. 95-00268 from the United States- 
Israel Binational Science Foundation (BSF), Jerusalem, Israel. D. Kandel is the incumbent 
of the Ruth Epstein Recu Career Development Chair. 

APPENDIX A: 

Here we study the linear stability of circular, uniformly spaced steps with unit step 
separation in the absence of step-step interactions. This configuration with p n = n is not a 
steady state. However, in the large radius limit the steps move very slowly, and the uniform 
state is extremely close to a steady state. We therefore refer to it as a quasi steady state. 
We regard the quasi steady state as unstable if the growth of perturbations is faster than 
the motion of steps. 
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In the absence of step-step interactions Eqs. (|10D simplify to 



p n = — = , with (Al) 

dr p n 



1 1 



Pn Pn + 1 



(l-g)ln^ -q{- + — 

y ^' Pn + 1 ^ \Pn Pn+1, 

(A2) 

For large values of n, the velocity of the nth step in the quasi steady state is given by: 



n 3 



Pn = — + O (n- 4 . (A3) 
l + q v 7 

To check the linear stability of the above configuration we perturb the step positions ac- 
cording to 

p n = n + Ae 1 ^-^ . (A4) 

Equating the time derivative of this perturbation with the velocity of the nth step results 
in an equation for to. We find that to the lowest orders in A and n" 1 

4ig(l - cos0)n~ 2 
{l + q) 

Since the magnitude of to is significantly larger than the step velocities in the quasi steady 
state, positive values of Im{uo) lead to an instability We see that away from the DL case 
(q > 0) the system is unstable. The most unstable mode is = —tt (step pairing). The 
= mode (uniform translation) is marginal as is the q = (DL) case. 



APPENDIX B: 

In this appendix we give some technical details of the algebraic manipulations performed 
in sections IV and V. 

To calculate expressions A and B in Eq. ([25]), we first express the scaled positions of the 
steps as power series in 9~ x . Inserting Eq. (|24J) into Eq. (^) we find that to fifth order in 

e- 1 

5 B~ l 3 F' 6~ 2 5 (-3 F' 2 + F F") 6~ 3 (-93 F' 3 + 61 F F' F" - 6 F 2 F&) 



Xi_2 = X — ; — 1 ^ — — h 



2F F 3 2F 5 AF 7 

1935 F' 4 - 1890 FF' 2 F" + 260 F 2 F' F^ + F 2 (l80 F" 2 - 17 F F^)) 6~ 5 

24 F 9 

3d -i F > 9 -2 ^F' 2 + FF")e- 3 (33 F' 3 -21FF' F" + 2F 2 F^ 0- 



+ 



2F F 3 2F 5 12 F 7 

135 F' 4 + 126 F F' 2 F" - 12 F 2 F" 2 - 16 F 2 F' F® + F 3 F&) 



24 F 9 

Xi X 2F 
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e- 1 

1+1 2F 

3 0-1 F , 0-2 (_3 F' 2 + F F") 6~ 3 (33 F' 3 — 21FF' F" + 2 F 2 F< 3 >) 6^ 

X i+2 = X + — — — Y^fl 

-135 F /4 + 126 F F' 2 F" - 12 F 2 F" 2 - 16 F 2 F' F® + F 3 F< 4 ') r 5 



— # + 



24 F 9 

5 0- 1 3 F' fl~ 2 5 (-3 F' 2 + F F") 6~ 3 (-93 F /3 + 61 F F F" — 6 F 2 F&) 6^ 



2F F 3 2F 5 4F 7 



f 1935 F /4 - 1890 FF' 2 F" + 260 F 2 F' F^ + F 2 ( 180 F" 2 - 17 F F< 4 >)) 

+ 5 " , (Bl) 

24 F 9 > \ J 

where F and its derivatives are evaluated at x — (x, + x i+ i) /2. 

Using these expressions we expand the step velocities in 9~ l and obtain Eq. (p5|). In the 
general case the expressions for A and B are too cumbersome to be written here. Instead 
we give these expressions in the two limiting cases. In the DL case 

(2F AF' 10 F 2 10 F" 18 F' F" •> 
^dl = 5 ( 5- + "TT— + + F + 6F® + 



i 3 f Fi x F J F x 3 

3 5F' 5F /3 7F" 10 F'F" 9 F" 2 F' 2 (5 F + 9 x 2 F") 



Bdl = 9 7 + "7; q - -=5 + —7^ + 



v x 4 Fx 3 F 3 x Fx 2 F 2 x F 2 F 3 x 2 

5j pO) QF t F {3) 3F (4)\ 3 F' 

+ _ f7~ + F 2 + ~F~ J ~ F 2 ^ 4 ~ F 3 ^ 3 ' ( 2) 

Inserting Eq. ( p2]) into Eq. (|28f) we find the differential equation which govern the scaling 
function in the DL case: 



g \12 F F {3) + 3 F F {4) + 9 F" 2 



15 F' F" + 5FF® 



x 

x z x 6 x 4 / x 4 4 

In the ADL case 

(I F' 3 F' 2 3 F' 3 5 F" 6 F' F" 3 F^ \ 1 F' 
Aadl -9 p—p + -pT^ ~ ~F 3 ~ + ~Fx~ + F 2 + _ F~ ) + F 2 ^ 3 + F^x 2 ' 

/ 3 F F' 2 3F' 3 9F /4 3F" F' F" 21 F' 2 F" 

&ADL — 9 — ~ ^ 7 + o — ^ n ^ 1" - ^ 9 + 



2Fx 4 F 2 x 3 F 3 x 2 F 4 x 2F 5 F 2 x 2 2F 3 x 2F 4 
3F" 2 5F( 3 ) 3F'F( 3 ) 3F^\ 3 2F' 3F 2 F" 

+ F 3 + 2F 2 x + 2F 3 + 2F 2 J ~ 2F 3 x 4 ~ F 4 x 3 ~ 2F 5 x 2 + 2F 4 x 2 ' ( ^ 

and the differential equation for the ADL scaling function is: 

/ 3F 3F' 3F' 2 3F /3 3F' 4 3F" 3 F'F" 15 F' 2 F" 3 F" 2 5F® 
9 \~2~x 1 + 2~^~ 2Fx 2 ~~ 2F 2 x + ~F 3 x 2 " + Fx 2F 2 + F + 2x 

3F / F (3) 3j p(4)\ 3 3F / xF > pn 

+ F + 2 J ~ 2Fx 4 ~ 2F 2 x 3 ~ ~ F 3 x 2 + 2F 2 x 2 ~° ' 
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In section V we used the moments 

Mo = f (^ + F 2 B ) jdx, 

M 2 = Jx 2 + F 2 Bj dx . 

to set boundary conditions for the scaling function at x . 
In the DL case 

/ F 2 2 F F' 5 F' 2 5 F F" \ 1 

Mo DL = 9 (^- 2 -^ + — + — + 9F>F" + 3FF>" 

\ X a X A X X J X 6 

M g (¥H _ 6F F' - xF' 2 - xF F" + 9x 2 F' F" + 3x 2 F F'") + - . 

V x x 



In the ADL case 



M 0aol = 9 - 77^ + 



F F' 3F' 2 3F' 3 5F" 3F' F" 3F 



2x 3 2x 2 2xF 2F 2 2x 
1 F' 



2x 3 F 2x 2 F 2 ' 
(3F 5F' 3xF' 2 3x 2 F' 3 xF" 3x 2 F' F" 3x 2 F" r 

2adl ~ 9 \2x^~^2 2F 2~F 2 2~ + F + ~^2~ t 

3 F' 

+ 2x~F + 2F 2 ' 



APPENDIX C: 

In this appendix we show that in the scaling state, the step density near the facet edge 
must vanish when the facet size diverges. This implies that F (x ) = where x is the scaled 
position of the facet edge. We restrict ourselves to situations where steps collapse towards 
the origin one at a time (consistently with simulation results). 

It is tempting to argue that since our system is expanding and slowing down, in the long 
time limit it approaches the equilibrium state of straight steps in contact with a facet. In 
this equilibrium state the step density near the facet edge vanishes as a square root. This 
equilibrium state is not consistent with our density function, which also predicts a square 
root approach to zero, but with a time dependent coefficient. Moreover even in the long 
time limit steps in our system continue to collapse. In each collapse period, just before the 
step disappears, its chemical potential diverges. At these times the system is not close to 
equilibrium. 

Consider the velocity of the first step. It depends on the positions of the first three steps 
Pi, p 2 and p 3 through 

l £i — £2 / 01 \ 

Pi = — ' 7 r • CI 

Pi (l-q)]n^-q(± + ±) 

p2 is always larger then p 1 so the denominator is always negative. The direction of motion 
of the first step is thus given by the sign of the numerator 
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* r = h-Z* = T-T +29 (r J -^- n+n ( 1 v ) > (C2) 
Pl P2 V (P2 - Pi) P2 + P3 (p 3 - P2) J 



which must be positive to avoid bounding of the first and second steps. 
Fixing p 2 and p 3 , the value of p* which minimizes Af is found by solving 



6g 1 



dpi 

The only solution between p 2 and the origin is 



(C3) 




This is indeed a minimum since the second derivative of M with respect to p\ is always 
positive when < p\ < p 2 . 

Substituting pl into M we find that 

M (pl) = „ \ - i + 



P2 + y¥-^%+f ^ (y'p^+f -yf) 3 (P2 + P3)(P3-p 2 ) 3 ' 

(C4) 

The first three terms of M {p\) decay when p 2 is large. The fourth term however will remain 
finite unless the difference p 3 — p 2 diverges with p 2 . If the fourth term does remain finite, 
■Af (p*) becomes negative, the velocity of the first step becomes positive and the first step 
cannot pass p* on its way to the origin. 

The above analysis indicates that the step density vanishes near the facet edge. We 
consider two scenarios for the collapse of the first step. In the first scenario the step starts 
the collapse from a position below p\. In this case the separation between the first and 
second steps is initially divergent since p 2 — p\ diverges as y/p~2. But the initial configuration 
in the collapse of one step is the final configuration of the former collapse period. Thus the 
former period ended with a divergent distance p 3 — p 2 . 

In the second scenario the first step starts to collapse from a position above p*. If the 
first step collapses alone (i.e p 2 is large throughout the collapse period), the distance p 3 — p 2 
must grow to allow pi to pass through p*. 

Either way there must be some point in the collapse period where the separation p 3 — p 2 
diverges. At this point the step density near the facet edge vanishes. 
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